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Abstract 

We consider general properties of groups of interacting oscillators, for which the natural fre- 
quencies are not in resonance. Such groups interact via non-oscillating collective variables like 
the amplitudes of the order parameters defined for each group. We treat the phase dynamics of 
the groups using the Ott-Antonsen ansatz and reduce it to a system of coupled equations for the 
order parameters. We describe different regimes of co-synchrony in the groups. For a large number 
of groups, heteroclinic cycles, corresponding to a sequental synchronous activity of groups, and 
chaotic states, where the order parameters oscillate irregularly, are possible. 

PACS numbers: 05.45.Xt, 05.45.Ac 
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I. INTRODUCTION 



Models of coupled limit cycle oscillators are widely used to describe self-synchronization 
phenomena in various branches of science. The applications include physical systems like 

n n n 

Josephson junctions lasers [2j, and electrochemical oscillators [3], but similar models 
are also used for neuronal ensembles [4|, the dynamics of pedestrians on bridges jl], 6], 
applauding persons [7], etc. 

In many cases a model of a fully connected (globally coupled) network is appropriate, it 
means that the oscillator population is treated in the mean field approximation. Ensembles 

sustained oscillators are successfully handled in the framework of 



of weakly interacting sel 

phase approximation js-ll2|. Most popular are the Kuramoto model of sine-coupled phase 
oscillators, and its extension, the Kuramoto-Sakaguchi model [131 ] . This model describes self- 
synchronization and appearance of a collective mode (mean field) in an ensemble of generally 
non-identical elements as a nonequilibrium phase transition. The basic assumptions behind 
the Kuramoto model are that of weak coupling and of closeness of frequencies of oscillators, 
the latter results in the presence of resonant terms in the coupling function only. References 
to detailed aspects of the Kuramoto model can be found in [14l4l6] . 

In many cases the ensembles of oscillators are not uniform and can be considered as 
consisting of several subensembles (e.g., in brain different groups of neurons can have differ- 
ent characteristic rhythms). If one still assumes that the frequencies of these subgroups are 



close (compared to the coupling), then a model of several interacting subpopulations 17. 18| 
or of an ensemble having a bimodal (or a multi-modal) distribution of frequencies 19M25| 
is adopted. Similarly, one can also model two ensembles, one of which consists of active 
and another of passive elements, which are coupled resonantly due to closeness of their 



frequencies |26|, 



27] 



In this paper we study a novel situation of non-res onantly coupled oscillator ensembles. 
We assume that there are several groups of oscillators, the frequencies in each group are close 
to each other, but are strongly (compared to the coupling strength) different between the 
groups. In this situation the coupling within the group is resonant, like in usual Kuramoto- 



type models, but the coupling between the groups can be only non-resonant 28]. It means 
that the coupling can be via non-oscillating, slow variables only, i.e. via the amplitudes 
of the mean fields. In the context of a single Kuramoto model such a dependence on the 



amplitude of the mean field corresponds to a nonlinearity of coupling, recently studied 



m 



29M32| . Nonlinearity in this context means that the effect of the collective mode on 



an individual unit depends on the amplitude of this mode, so that, e.g., the interaction of 
the field and of a unit can be attractive for a weak field and repulsive for a strong one. 
Mathematically, this is represented by the dependence of the parameters of the Kuramoto- 
Sakaguchi model (the coupling strength, the effective frequency spreading, and the phase 
shift) on the mean field amplitude. Here we generalize this approach to several ensembles, so 
that the parameters of the Kuramoto-Sakaguchi model describing each subgroup depend on 
the mean field amplitudes of other subgroups (e.g., resonant interactions within a group of 
oscillators can be attractive or repulsive dependent on the amplitude of the order parameter 
of another group). 

In Section [Til we introduce the basic model of non-resonantly interacting ensembles. We 
also formulate the equations for the mean fields of the ensembles following the Ott-Antonsen 



theory |33j,|34|. The simplest situation of two interacting ensembles is studied in Section [TTTl 



In Section |IV] we describe three and several interacting ensembles, focusing on nontrivial 
regimes of sequential synchronous activity following a heteroclinic cycle, and on chaotic 
dynamics. 



II. BASIC MODEL OF NON-RESONANTLY INTERACTING OSCILLATOR EN- 
SEMBLES 

A. Kuramoto-Sakaguchi model and Ott-Antonsen equations for its dynamics 

A popular model describing resonant interactions in an ensemble of oscillators having 



close frequencies is due to Kuramoto and Sakaguchi 13] 



j> k = u k + lm{KZe~ i * k ), Z = J^H^^ k = l,...,N. (1) 

Here <pk is oscillator's phase, Z is the complex order parameter (mean field) that also serves 
as a measure for synchrony in the ensemble, u k are natural frequencies of oscillators, and 
K = 2a + 2ib is a complex coupling constant. Recently, Ott and Antonsen 33j, |34| have 



demonstrated that in the thermodynamic limit N — > oo, and asymptotically for large times 
the evolution of the order parameter Z in the case of a Lorentzian distribution of natural 



frequencies g(co) = A[tt(uj — ujq) 2 + A 2 ] 1 around the central frequency ujq is governed by a 
simple ordinary differential equation 



Written for the amplitude and the phase of the order parameter defined according to Z = 
pe 1 ®, the Ott-Antonsen equations 



are easy to study: Eq. ([3]) defines the stationary amplitude of the mean field (which is non- 
zero above the synchronization threshold a c = A), while Eq. (j3J) yields the frequency of the 
mean field. 

B. Non-resonantly interacting ensembles 

We consider several ensembles of oscillators, each characterized by its own parameters 
co , A, a, b. The main assumption is that the central frequencies u of different populations 
are not close to each other, and also high-order resonances between them are not present. 
Such a situation appears typical for neural ensembles, where different areas of brain demon- 
strate oscillations in a very broad range of frequencies, from alpha to gamma rhythms. 
Because there is no resonant interaction between the oscillators in different ensembles, they 
can interact only non-resonantly, via the absolute values of the mean fields. Assuming 
that only Kuramoto order parameters ([1]) (but not higher-order Daido order parameters 
Z m = (e m ^)) enter the coupling, a general non-resonant interaction between populations 
can be described by the dependencies of the parameters ojq, A, a, b on the amplitudes of the 
mean fields pi, where index I counts the subpopulations. Moreover, one can see from (I3f4"l) 
that the equation for the amplitude is independent on the phase, therefore we can restrict 
our attention to the amplitude dynamics ([3]). Furthermore, we assume the coupling to be 
week, so only the leading order corrections ~ p 2 are included. All this leads to the following 
general model for interacting populations 




(2) 



P 



-Ap + o(l- p 2 )p , 
lu + 2bp 2 , 



(3) 
(4) 



pi = (-Ai - T lm p 2 m )pi + (a t + A lm p 2 m ){\ - p 2 )p h 



/ = !,.. .,L 



(5) 
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with coupling constants Ti m , A\ m . Note that because the widths of the frequencies distribu- 
tion cannot be negative, coefficients Ti m must satisfy A^ + Ti m > 0. Below we assume that 
there is no nonlinearity inside ensembles Tu = An = 0. 

In this paper we will not investigate model (J5J) in its full generality, as it would require 
a rather tedious analysis. Instead, we will consider two simpler models, which describe par- 
ticular types of interaction, but nevertheless allow us to demonstrate interesting dynamical 
patterns. In model A we assume that only frequencies are influenced by the coupling, i.e. 
Ai m = 0. This leads to a system 

pi = (ai - Ai - Ti m p m - aipf)pi (6) 

Another model B takes into account the interaction via coupling constants only (i.e. Ti m = 
0); additionally we will assume here that the distributions of frequencies in all interacting 
ensembles are narrow A; — > 0. In the limit of identical oscillators we obtain from flSJ) 

pi = {ai + A lm p 2 m ){l -pf)pi (7) 

Here we note that the Ott-Antonsen equations for the ensemble of identical oscillators de- 
scribe not a general case, but a particular solution, while a general description delivers the 



Watanabe-Strogatz theory [35|, |36|. Thus the dynamics of model B should be considered as 



a special singular limit A — > 0. 

Below, in sections DTD and [TV] we describe the dynamics of these two models, for the cases 
of two, and three and more interacting ensembles, respectively. 



III. TWO INTERACTING ENSEMBLES 

Let us first rewrite models (I61PTI) for the simplest case of only two interacting ensem- 
bles. Additionally, for model A we assume ai = 1 (equivalently, one could renormalize the 
amplitudes of order parameters p\^ to get rid of these coefficients). Thus the model A reads 

Pi = Pi(5i -d l2 p 2 2 - p\) , 
P2 = P2($2 ~ d 2 ip\ - pi) ■ 
For model B a normalization of amplitudes is not possible, and it reads 

Pl = e lPl {l - D 12 p 2 2 ){l - pi) , 

(9) 

P2 = e 2 pi(l - D 21 p\)(l - p\) . 
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Generally, parameters 5i )2 = 1 — ^1,2, dik = T^, ei )2 = ai,2, Dik = —■^%kl a i can have 
different signs. 

As the first property of both models we mention that the dynamics is restricted to the 
domain < pi j2 < 1. Formally, this follows directly from (J5J), physically this corresponds to 
the admissible range of values of the order parameter. Furthermore, for model A (jSJ) we can 
apply the Bendixon-Dulac criterion 

» I 1 .A + JL (_LA = ^1±A < 



dpi \p1P2 J dp 2 \P1P2 J P1P2 
from which it follows that it cannot possess periodic orbits. 

Remarkably, model B (jHJ) can be written as a Hamiltonian one. With an ansatz 

expt/1,2 =p?,2(l-Pi,2) -1 ( 10 ) 

it can be represented in a Hamiltonian form 

dH(y 1 ,y 2 ) . dH{y 1 ,y 2 ) 
Vi = « , V2 = ~ , 

H = 2e iy2 - 2e 2Vl - 2e x D X2 ln(l + e^ 2 ) + 2e 2 D 21 ln(l + e yi ) . 
Thus model B may demonstrate a family of periodic orbits if the levels of the Hamiltonian 
are closed curves. We stress that the Hamiltonian structure of the model does not exclude 
existence of stable equilibria at p = 0, 1 because the transformation (TTOl) is singular at these 
states; in the Hamiltonian formulation ( ITT]) these stable equilibria correspond to trajectories 
moving toward =Foo. 

The dynamics of both models is mainly determined by the existence and stability of 
equilibria. For model A (jHJ) possible equilibria are the trivial one Si (0,0), two states where 

1/2 1 /2 

one of the order parameters vanish S 2 (5 1 , 0) and S^O, 5 2 ), and a state where both order 
parameters are non-zero ^((^i — di 2 5\){l — di 2 d 2 i)^, (5 2 — d 2 iSi)(l — dud^)^ 1 ). Similarly, 
model B © always has equilibria Mi(0, 0), M 2 (l, 0), M 3 (0, 1) and M 4 (l, 1), and additionally 
a nontrivial state M 5 (D^" 1 1 ^ 2 , D^^ 2 ) existing if Z) 12 , D 2 i > 1. 

We illustrate possible types of dynamics (up to symmetry 1 H 2) in models A,B in 
Figs. [TH21 Here it is worth mentioning, that model (jHJ) is structurally of the same type 
as typical models of interacting populations in mathematical ecology {3?]]. Model B (J9j) 
resembles them as well, but has a distinctive property that fully synchronized cluster p — 1 
is invariant. Referring for the details to Appendix |Aj we describe briefly possible regimes in 
these models. 
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FIG. 1. (Color online) Six different patterns of the dynamics of system ([8]). (a) Global stability 
of a trivial state (for Si 2 < 0). (b,c) Global stability of 54 when both populations are partially 
synchronous (conditions for this are (|A2p for (b) or ()A3|) for (c)). (d) Competition between clusters 
if the coupling is strongly suppressive (|A6P ); here we have bistablity of states 62,3 describing 
synchronous activity of one cluster and asynchronous of another one. (e) Asymmetric interaction 
between clusters arises under condition (|A8j) ; here always a heteroclinic trajectory from saddle 
point S3 to stable node S2 exists (red dashed line), (f) Global stability of S2 under condition (|A8|) , 

1. Global stability of trivial equilibrium point Si(0, 0), Mi(0, 0) (FigJT^i, Figf2^,b) means 
that a fully asynchronous state is stable in both ensembles. 

2. Stability of a nontrivial state off coordinate axes S4 and M 4 (FigJTJ^c, Figf2fc,d). Here 
both ensembles are synchronized (in model A not completely because of a distribu- 
tion of frequencies, in model B completely because we assume identical oscillators in 
ensembles) . 

3. Competition between ensembles (FigJT]ij2k): Only one ensemble synchronizes while 
the other one desynchronizes. Which ensemble is synchronous depends on initial con- 
ditions. 

4. Suppression: One ensemble always "wins" and is synchronous while the other one 
desynchronizes (steady states S2,M 2 are global attractors, of course also stability of 
"symmetric" states S 3 ,M 3 is possible) (FigJT^,f j2]F,g,h). 
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5. The case of bistability of the trivial and the fully synchronous states of both ensembles 
(Figj2j) is possible in the model B only. 

6. Periodic behavior (FigJ2J) is possible only in ensemble B, it corresponds to an inter- 
action of populations of "predator-pray" type. Because of the system is Hamiltonian, 
the oscillations are conservative like in the Lottka-Volterra system. 

While in our analysis we studied models (!H|9|) describing dynamics of the order parameters 
in the Ott-Antonsen ansatz, all the regimes described above can be observed when one 
simulates original equations of the ensembles of phase oscillators 01) , at sufficiently large 
number of units N. In Fig. [3] we illustrate two nontrivial regimes of two subpopulations 
of phase oscillators at N = 10 3 . Figure |3ja) shows the dynamics of mean fields in the 
case of a competition between two subpopulations that interact via frequency mismatch 
modulation, see Fig. [TJd). Figure [3](b) illustrates a periodic behavior of two subpopulations 
like in Fig.Eft). 

IV. THREE AND MORE INTERACTING ENSEMBLES 

In this section we generalize the results of Section II I II to many interacting ensembles. 
We do not aim here at the full generality, but rather present interesting regimes based on 
the elementary dynamics depicted in Figs. [Tf2l According to the consideration above, we 
restrict our attention to two basic models A and B (j7|). Generally, model B cannot be 
rewritten in a Hamiltonian form, but by applying transformation (jiup one can easily see 
that this system has a Liouvillian property - the phase volume is conserved. 

A. Symmetric case: cosynchrony and competition 

Here we describe mostly simple regimes that are observed in a symmetric case where 
parameters of all ensembles and their interaction are equal. This corresponds to equal values 
a>i = a, Ai = A, Ti m = V in and a\ = a, A[ m = A in (j7|). In model A, the only nontrivial 
regimes are those where asynchronous states are unstable A < a. Then one observes either 
a coexistence of synchrony like in Fig. [Tb (for T < a) or a competition like in Fig. [TH (for 
r > a). In the latter case only one ensemble is synchronous, while other desynchronize. 
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FIG. 2. (Color online) Ten different dynamical regimes in system Q. (a,b): Global stability of 
the trivial state, arises at conditions (|A1|) . Case (a): D12 < 1, case (b): D\i > 1. (c,d): Global 
stability of M±(l, 1) when both clusters are in the synchronized state, under condition of a weak 
suppressive coupling (|A4|) for (c) or at (|A5|) for (d). (e): Competition between clusters, arises 
at strong suppressive coupling (|A7[) . Here we have bistablity of steady states M2 5 3; each of these 
points corresponds to synchronous activity of one cluster and asynchronous of another one. (f): 
Asymmetric interaction between clusters at asymmetric coupling (lAlljl . Here always a sequence 
of heteroclinic trajectories M3 — > M4 — > M2 (red dashed lines) is present. (g,h): The situation 
of global stability of fixed point M3 while conditions ()A10|) are satisfied (case (g): D%\ > 1, case 
(h): D21 < 1). (i): Bistablity of fully asynchronous and fully synchronous states, arises if (|A12p 
is valid. In this case stable manifolds of the saddle point M5 divide basins of attraction of stable 



points Mi, M4. (j): The case of periodic behavior, arises at conditions (|A13j) , 

9 



1 

0.8 
0.6 
0.4 
0.2 
0, 

















... Pl 




" 1 

1 
t 








— p 2 




1 










/ 











200 400 600 800 1000 
time 



0.6 
0.4 
0.2 




J \L 



ntj r 



(b) 



2000 4000 6000 8000 10000 
time 



FIG. 3. (Color online) Modeling of ensemble consisting of two subpopulations of N = 10 3 phase 
oscillators, (a) Subpopulations interact via modulation of effective frequency mismatch {SJ. Case 
of competition between subpopulations for parameter values #1,2 = 10, d\2 = ^21 = 12. (b) 
Subpopulations interact via coupling modulations (|9|). A periodic regime is presented at parameter 
values Ei = —1, £2 = 1, D\2 = D21 = 2. To avoid a spurious clustering and to ensure validity 
of Ott-Antonsen description, a small mismatch was added: co n were randomly distributed in the 
range [-0.025, +0.025]. 




(b) 




FIG. 4. (Color online) Multistability of steady states C n corresponding to synchronous state of 
only one cluster for (a) system (jSJ) (A > a, T > a) and (b) system ([7]) (a > 0, A < —a). 



Similar regimes can be observed in model B for a > 0, A < —a. Additionally, in model B 
a coexistence of full synchrony in all ensembles and a full asynchrony, like in Fig. |2} can be 
observed for a < 0, A > —j^j- We illustrate the regimes of competition in Fig. H]for the 
case of three interacting populations. 
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B. Heteroclinic synchrony cycle 



Here we discuss a multidimensional generalization of the interactions where in a pair of en- 
sembles one group always synchronizes while another one is asynchronous (see Figs.[T](e)J2(f)). 
In the examples presented in these graphs, both ensembles would self-synchronize separately, 
but due to interaction synchrony in ensemble 2 disappears while ensemble 1 remains syn- 
chronous. One can say that in synchrony competition between the first and the second 
ensembles, the first ensemble wins. Suppose now, that a third self-synchronizing ensemble 
is added, which wins in the competition with the first one but looses in the competition to 
the second one. Then a cycle 2— > 1 — ^ 3 — > 2 — >• 1 — > 3 . . . will be observed. Moreover, 
because in the dynamics Figs. (He)j2](f) the transition 2 — > 1 follows the heteroclinic orbit 
connecting steady states S3 and S2, the cycle in the system of three ensembles will be 
a heteroclinic one, with asymptotically infinite period. Such a cycle has been studied in 
different contexts 38N40| . For a review of robust heteroclinic cycles see 4ll. l42lj (sometimes 



one uses a term "winnerless competition" to describe such a dynamics [43 



We demonstrate the heteroclinic synchrony cycle for three interacting ensembles in Fig. [51 
One can see that synchronous states of ensembles appear for longer and longer time inter- 
vals. It is interesting to note that heteroclinic cycles have been observed in ensembles of 



identical coupled oscillators 45M50I]. There the nontrivial dynamics is in the switchings of 
full synchrony between different clusters. In this respect the heteroclinic cycle in the model 
B resembles such a regime. On the other hand, the heteroclinic cycle in model A is different: 
here the natural frequencies of oscillators are different and the states of synchrony are not 
complete, so the identical clusters never appear. 

Finite size effects are nontrivial for the heteroclinic cycles described. Indeed, it is known 
that while in the thermodynamic limit deterministic equations for the order parameters 
can be used, finite size effects can be modeled via noisy terms that scale roughly as ~ 



N- 1 / 2 [51H53J. On the other hand, noisy terms destroy perfect heteroclinic orbit, making 
the transitions times between the states finite and irregular. Exactly this is observed at 
modeling the interacting finite size ensembles (Fig. [6]). While for small N the heteroclinic 
cycle is completely destroyed, for large N it looks like a noisy limit cycle. 
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FIG. 5. (Color online) Stable heteroclinic cycles caused by asymmetric interactions between clusters 
in system ([6]) (a,c) and in system (JTj) (b,d). Parameters: (a,c) ai — 6i > 0, Ti2 > 

r« > ^fer 1 ' r ^3 > ^fe?, r 21 < s^l, r 13 < ^r 1 * r 3 2 < ^te? 1 ^ (d) > o, 

^4i2 < —a%, A 31 < -a 3 , A 23 < -a 2 , A 2 \ > -a 2 , Ai 3 > -ax, A 32 > -a 3 . Panels (a,b) show the 
phase space portraits while time series are presented in panels (c,d). 




FIG. 6. (Color online) Dynamics of the order parameters of three interacting populations of 
oscillators (parameters like in Fig. [5] (a,c)) for three different sizes of populations: (a) N = 100, 
(b) N = 400, and (c) N = 10000. 
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C. Chaotic oscillations 



Here we discuss possible "predator-pray" -type regimes (cf. Fig. [2]) for many ensembles. 
An elementary "oscillator" depicted in Fig. |2j can be represented as a Hamiltonian system 
with one degree of freedom. Several of such elementary conservative "oscillators", being 
coupled, can yield quasiperiodic and chaotic regimes. In the case of two interacting con- 
servative "oscillators" (i.e. of four interacting ensembles), system (j7|) can be rewritten as 
follows: 

Pl,2 = ei, 2 pi, 2 (l - A)P2,1 - £>H>?, 2 )(1 - Pl,2) , , , 

(12) 

V\,2 = £1,2^1,2(1 - A)^2,l - ^ ) lPl )2 )( 1 - V l,2) ■ 

Here the parameters of the system were chosen in such a way that each pair of subpopulation 
(pi, p 2 ) an d (i>i, d 2 ) exhibits periodic oscillation being decoupled from another pair (at D\ = 
0), i.e. E\E2 < and D > 1. When the coupling between the two pairs is introduced 
(i.e. D\ 7^ 0), then in dependence on this coupling and initial conditions the dynamics 
can be qusiperiodic or chaotic. Like in general Hamiltonian systems with two degrees of 
freedom, it is convenient to represent the dynamics as a two-dimensional Poincare map. As 
a Poincare section (Figure [7^) we have taken the plane (vi, i; 2 ) at moments of time at which 
the variable p\(t) has a maximum. At small values of the coupling between the "oscillators" 
D\ the dynamics is typically quasiperiodic. While increasing Di, one can observe a transition 
to dominance of chaotic regimes in the system (Tl2|) (see FigJ7^,b and calculation of Lyapunov 
exponents in FigJTfc). Furthermore, we have confirmed the existence of chaotic oscillations 
by direct numerical simulation of four subpopulations satisfying f fl2|) . consisting of N = 10 3 
elements each (FigJTU). 



V. CONCLUSION 

In this paper we have introduced and studied a model of non-resonantly coupled ensembles 
of oscillators. It is assumed that oscillators form several groups, in each group the natural 
frequencies are close to each other, but the frequencies of different groups are rather different. 
This means that only oscillators within each group interact resonantly (i.e. the coupling 
terms depend on their phases), while interactions between the groups can be only non- 
resonant, i.e. depending on slow non-oscillating variables only. As a particular realization 
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FIG. 7. (a) Poincare sections on the plane (t>i,f2) demonstarting regular and chaotic dynamics 
at different values of D\ in the system (|12p . (b) Time series of a chaotic regime of system (|12p . 
for parameter values D\ = 0.5, Dq = 2.0, e\ = —1.0, £2 = 1.0. (c) Lyapunov exponents calculated 
at different values of D\, for some particular value of the Hamiltonian. From four Lyapunov 
exponents two always vanish, while other two vanish for small D\ (quasiperiodicity) and are non- 
zero for larger couling (chaos), (d) Chaotic time series of order parameters of four subpopulations 
of oscillators consisting of N = 10 3 elements each (the coupling configuration and the parameters 
are like in panel (b)). 

of such a setup we considered phase oscillators, which resonantly interact according to 
the Kuramoto-Sakaguchi model, and the non-resonant terms appear as dependencies of 
the parameters of the Kuramoto-Sakaguchi model on the amplitudes of the mean fields 
(Kuramoto order parameters) of other groups. 

We employed the Ott-Antonsen theory allowing us to write a closed system of equation 
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for the amplitudes of the order parameters. Analysis of this system constitutes the main 
part of the paper. The system resembles the Lottka-Volterra type equations used in math- 
ematical ecology for the dynamics of populations, but has nevertheless some peculiarities. 
For two coupled ensembles we demonstrated a variety of possible regimes: coexistence and 
bistability of synchronous states, as well as periodic oscillations. For a larger number of 
interacting groups more complex states appear: a stable heteroclinic cycle and a chaotic 
regime. Heteroclinic cycle means a sequence of synchronous epochs that become longer and 
longer. In a chaotic regime the order parameters demonstrate low- dimensional chaos. While 
the main analysis is performed for the Ott-Antonsen equations that are valid in the ther- 
modynamic limit of infinite number of oscillators in ensembles, we have checked finite-size 
effects in several regimes by modeling finite ensembles. Finiteness of ensembles only slightly 
influences the dynamics in most of the observed states, except for the heteroclinic cycle. 
Here a small effective noise due to finite-size effects destroys the cycle, producing nearly 
periodic noise-induced oscillations. 

One of the models we studied was that of groups of identical oscillators. Here in many 
cases only the states where some groups completely synchronize (i.e. all oscillators form an 
identical cluster) while other completely desynchronize (order parameter vanish) are possible. 
Heteroclinic cycle in this model also connects such states. There is, however, a nontrivial set 
of parameters, at which the order parameters of ensembles oscillate between zero and one, 
thus demonstrating time-dependent partial synchronization. Moreover, for four ensembles 
these oscillations are chaotic. This regime is quite interesting for a general theory of collective 
chaos in oscillator populations (cf. chaotic dynamics of the order parameter in an ensemble 



of Josephson junctions reported in [35j) and certainly deserves further investigation. 
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Appendix A: Details of analysis of two interacting ensembles 

Here we present details of the analysis of models (jHED, giving the conditions for different 
regimes presented in Figs. [Q [2j 

1. The case of global stability of trivial equilibrium point 5*1(0,0), Mi (0,0) (Figfjji, 
Figfjk^b). For system (JHJ) such a situation occurs in the case 5\ t2 < 0. For system (J9]) 
global stability of trivial state M x occurs if e^ 2 are negative and at least one of D 12 or 
D 21 less than 1: 

£i, 2 <0 ) mm(D 12 , D 21 ) < 1 (Al) 

2. The case of stability of non-trivial state off coordinate axes S4 and M4 (FigJT|D,c, 
Fig]2b,d). For system (jSJ) this situation occur in two cases. The first situation appears 
if Sx t 2 > (when isolated subpopulations tends to synchrony) and suppressive couplings 
are weak: 

Ai = j5 2 - M21 > 0, A 2 = - 5 2 d 12 > 0. (A2) 

The states S 2s 3 have eigenvalues —5i5\, Ai, and —S 2 5 2 , A2, respectively, and therefore 
are saddles. The origin is an unstable node (5i i2 > 0) and therefore the state 5*4 is an 
attractor (note that S4 always exists while ( ]A2j) holds). We call this situation "case 
of weak suppressive couplings" because ( 1A2I) can be written as d\ 2 < d 2 i < y|j-. 

The second situation appears if one of the subpopulations approaches to the asyn- 
chronous state (negative 5) while another group tends to synchrony and has positive 
influence on the first subpopulation: 

<5i > 0, 5 2 < 0, Ai > or 5 l < 0, 5 2 > 0, A 2 > 0. (A3) 

Condition Ai > is equivalent to d 2 \ < what means that coupling d 2 \ should 
be negative and absolute value of d 2 \ should be large enough to maintain partially 
synchronous state inside the second cluster (positive influence). 

For the system (jUJ) the situation of global stability of M 4 can be produced by two types 
of conditions. The first case is that of positive e^ 2 and weak suppressive couplings: 

si,2 > 0, £i 2 < 1, D 21 < 1. (A4) 
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Another case of global stability of M4 occurs if 

£1 < 0, £ 2 > 0, D 12 > 1, D 21 < 1 or e 1 > 0, e 2 < 0, D 12 < 1, D 21 > 1. (A5) 

The latter case differs from the previous one only by the direction of the flow on lines 
Pi, 2 = and the type of unstable points Mi, M 2 , M 3 (FigfSH). 

3. The case of competition between subpopulations (FigJTHfSfe). In model (jSJ) this type 
of behavior arises when 

5i, 2 > 0, Ai < 0, A 2 < 0. (A6) 

According to flA6j) the points S 2 ,3 are stable, while Si is unstable node and S 4 is a 
saddle. This case corresponds to the situation of strong suppressive couplings 

5x8-1 5 2 5 2 

Competitive behavior in the system (jUJ) is produced by positive £1,2 and strong sup- 
pressive couplings between subpopulations: 

£0 > 0, D X2 > 1, D 21 > 1. (A7) 

4. The case of global stability of synchronous state of only one cluster (S 2 ,M 2 ) (FigJTfe,fj2]f,g,h). 

In model (jHJ) only one group is synchronous in two cases. The first trivial situation is 
similar to conditions ( 1A3j) (when one group approaches to asynchronous state while 
another one tends to synchrony) but in this case the active group does not have suffi- 
cient positive coupling to maintain synchronization in the asynchronous subpopulation 
(FigJTf): 

8x > 5 2 < Ai < or 8x < 5 2 > A 2 < 0. (A8) 



Under conditions ( )A8j) only one of the fixed points S2 or S 3 exists and Si is always 
unstable. The second case occurs at an asymmetric interaction of intrinsically active 
clusters (isolated clusters tend to synchronous regime): 

8 1)2 > 0, and Ai < A 2 > or Ai > A 2 < (A9) 



In other words, it appears when one coupling coefficient is strong enough to fully 

§2 S2 

2 Si ■ 
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suppress the synchrony in the opponent, for example d 2 i > y! 2 ) while another one is 



weak or even non-suppressing rf 12 < y|^. In this case one can prove that S 4 does not 
exist, point S2 is a stable node, S3 and Si are saddles. Thus all trajectories approach 
stable node S 2 which corresponds to the synchronous state of the first group and to 
the asynchronous state of the second one. Because of this on the plane (pi, p 2 ) always 
exists heteroclinic trajectory connecting saddle point S3 and stable equilibrium S2 (red 
line in FigfjJ)- 

Global stability of point M 2 (M 3 ) of system fl9]) occurs in several different cases. The 
first case is similar to the situation in the system ([8]) at conditions (IA8[) when one group 
tends to synchrony (5 n > 0), another one approaches trivial state (S m < 0) and syn- 
chronous group does not have sufficient positive infuence to maintain synchronization 
in the asynchronous group: 

£1 < 0, e 2 > 0, D 12 < 1 or e 1 > 0, £ 2 < 0, D 2 i < 1. (A10) 

Corresponding phase planes are presented in Figj2g,h. Another case is that of positive 
£1,2 > and asymmetric couplings: 

£1,2 > 0, D 12 > 1 D 2i < 1 or L>i2 < 1 D 2 i > 1. (All) 

Under described above conditions (lAlOj) . (IA1 1|) it is easy to show that only one stable 
fixed point M 2 (1,0) exists, so all trajectories approach M 2 . In the case of (lAllj) a 



sequence of heteroclinic orbits connecting M 2 and M 3 (red lines in Figj2]F) appears. 

5. The case of bistability of trivial and fully synchronous states (Figj2}). 

In model (Q this happens for negative £i j2 and strong synchronizing couplings: 

£1,2 < 0, D 12 > 1, D21 > 1. (A12) 

6. Periodic behavior (FigfJJ). 

In model ([9]) periodic solutions can be observed. Conditions 

D 12 > 1, £21 > 1, £ X E2 < (A13) 

provide saddle type of points Mi_4 and existence of equilibrium M5 with imaginary 
eigenvalues iiJ ^ ^ dl ld^ld^ 21 ~~^ • Because model ([9]) can be rewritten as a Hamiltonian 
one, one has a family of periodic orbits. 



18 



[1] K. Wiesenfeld and J. W. Swift, Phys. Rev. E 51, 1020 (1995). 

[2] A. F. Glova, Quantum Electronics 33, 283 (2003). 

[3] I. Kiss, Y. Zhai, and J. Hudson, Science 296, 1676 (2002). 

[4] D. Golomb, D. Hansel, and G. Mato, in N euro-informatics and Neural Modeling, edited by 
F. Moss and S. Gielen (Elsevier, Amsterdam, 2001), vol. 4 of Handbook of Biological Physics, 
pp. 887-968. 

[5] S. H. Strogatz, D. M. Abrams, A. McRobie, B. Eckhardt, and E. Ott, Nature 438, 43 (2005). 
[6] B. Eckhardt, E. Ott, S. H. Strogatz, D. M. Abrams, and A. McRobie, Phys. Rev. E 75, 021110 
(2007). 

[7] Z. Neda, E. Ravasz, Y. Brechet, T. Vicsek, and A.-L. Barabasi, Nature 403, 849 (2000). 

[8] Y. Kuramoto, in International Symposium on Mathematical Problems in Theoretical Physics, 
edited by H. Araki (Springer Lecture Notes Phys., v. 39, New York, 1975), p. 420. 

[9] Y. Kuramoto, Chemical Oscillations, Waves and Turbulence (Springer, Berlin, 1984). 
[10] H. Daido, Prog. Theor. Phys. 88, 1213 (1992). 
[11] H. Daido, Prog. Theor. Phys. 89, 929 (1993). 
[12] H. Daido, Physica D 91, 24 (1996). 

[13] H. Sakaguchi and Y. Kuramoto, Prog. Theor. Phys. 76, 576 (1986). 

[14] A. Pikovsky, M. Rosenblum, and J. Kurths, Synchronization. A Universal Concept in Non- 
linear Sciences. (Cambridge University Press, Cambridge, 2001). 

[15] J. A. Acebron, L. L. Bonilla, C. J. P. Vicente, F. Ritort, and R. Spigler, Rev. Mod. Phys. 77, 
137 (2005). 

[16] S. H. Strogatz, Physica D 143, 1 (2000). 

[17] N. Tukhlina and M. Rosenblum, J. Biol. Phys. 34, 301 (2008). 
[18] O. V. Popovych and P. A. Tass, Phys. Rev. E 82, 026204 (2010). 
[19] L. L. Bonilla, J. C. Neu, and R. Spigler, J. Stat. Phys. 67, 313 (1992). 
[20] J. D. Crawford, J. Stat. Phys. 74, 1047 (1994). 

[21] L. L. Bonilla, C. J. P. Vicente, and R. Spigler, Physica D: Nonlinear Phenomena 113, 79 
(1998). 

[22] L. L. Bonilla, Phys. Rev. E 62, 4862 (2000). 

19 



[23] E. Montbrio, D. Pazo, and J. Schmidt, Phys. Rev. E 74, 056201 (2006). 

[24] E. A. Martens, E. Barreto, S. H. Strogatz, E. Ott, P. So, and T. M. Antonsen, Phys. Rev. E 

79, 026204 (pages 11) (2009). 

[25] D. Pazo and E. Montbrio, Phys. Rev. E 80, 046215 (2009). 

[26] H. Daido and K. Nakanishi, Phys. Rev. Lett. 93, 104101 (2004). 

[27] D. Pazo and E. Montbrio, Phys. Rev. E 73, 055202 (2006). 

[28] Another novel type of interaction appears if the frequencies of two groups are in a high-order 



resonance like 2:1, see 



54]- 



[29] M. Rosenblum and A. Pikovsky, Phys. Rev. Lett. 98, 064101 (2007). 
[30] A. Pikovsky and M. Rosenblum, Physica D 238(1), 27 (2009). 
[31] G. Filatrella, N. F. Pedersen, and K. Wiesenfeld, Phys. Rev. E 75, 017201 (2007). 
[32] F. Giannuzzi, D. Marinazzo, G. Nardulli, M. Pellicoro, and S. Stramaglia, Phys. Rev. E 75, 
051104 (2007). 

[33] E. Ott and T. M. Antonsen, CHAOS 18, 037113 (2008). 

[34] E. Ott and T. M. Antonsen, CHAOS 19, 023117 (2009). 

[35] S. Watanabe and S. H. Strogatz, Physica D 74, 197 (1994). 

[36] A. Pikovsky and M. Rosenblum, Phys. Rev. Lett. 101, 264103 (2008). 

[37] J. D. Murray, Mathematical Biology. I. An Introduction (Springer, Berlin, 2002). 

[38] F. H. Busse and R. M. Clever, in Recent Development in Theoretical and Experimental Fluid 

Mechanics, edited by U. Mtiller, K. G. Roessner, and B. Schmidt (Springer, NY, 1979), pp. 

376-385. 

[39] T. Clune and E. Knobloch, Physica D 74, 151 (1994). 

[40] J. Guckenheimer and P. Holmes, Math. Proc. Camb. Phil. Soc. 103, 189 (1988). 
[41] M. Krupa, J. Nonlinear Sci. 7, 129 (1997). 

[42] V. Afraimovich, P. Ashwin, and V. Kirk, Editors, A focus issue on Robust Heteroclinic and 

Switching Dynamics, Dynamical Systems, vol. 25, n. 3 (2010). 
[43] M. Rabinovich, A. Volkovskii, P. Lecanda, R. Huerta, H. D. I. Abarbanel, and G. Laurent, 

Phys. Rev. Lett. 87, 068102 (2001). 
[44] V. Afraimovitch, I. Tristan, R. Huerta, and M. Rabinovich, CHAOS 18, 043103 (2008). 
[45] P. Ashwin and J. W. Swift, J. Nonlinear Sci. 2, 69 (1992). 
[46] D. Hansel, G. Mato, and C. Meunier, Phys. Rev. E. 48, 3470 (1993). 

20 



[47] H. Kori and Y. Kuramoto, Phys. Rev. E 63, 046214 (2001). 
[48] H. Kori, Phys. Rev. E 68, 021919 (2003). 
[49] P. Ashwin and J. Borresen, Phys. Rev. E 70, 026203 (2004). 
[50] P. Ashwin and J. Borresen, Physics Letters A 347, 208 (2005). 
[51] A. Pikovsky and S. Ruffo, Phys. Rev. E 59, 1633 (1999). 

[52] A. Pikovsky, A. Zaikin, and M. A. de la Casa, Phys. Rev. Lett. 88, 050601 (2002). 
[53] E. J. Hildebrand, M. A. Buice, and C. C. Chow, Phys. Rev. Lett. 98, 054101 (2007). 
[54] S. Luck and A. Pikovsky, Resonancly interacting multifrequency oscillator enesembles, in 
preparation (2011). 



21 



